﻿/*
 *  Copyright 2014-2022 The GmSSL Project. All Rights Reserved.
 *
 *  Licensed under the Apache License, Version 2.0 (the License); you may
 *  not use this file except in compliance with the License.
 *
 *  http:www.apache.org/licenses/LICENSE-2.0
 */

#include <stdio.h>
#include <string.h>
#include <assert.h>
#include <stdlib.h>
#include "sm2.h"
#include "sm2_common.h"



const peace_sm2_BN peace_sm2_P = {
	0xffffffffu,
	0xffffffffu,
	0x00000000u,
	0xffffffffu,
	0xffffffffu,
	0xffffffffu,
	0xffffffffu,
	0xfffffffeu,
};

const peace_sm2_BN peace_sm2_B = {
	0x4d940e93u,
	0xddbcbd41u,
	0x15ab8f92u,
	0xf39789f5u,
	0xcf6509a7u,
	0x4d5a9e4bu,
	0x9d9f5e34u,
	0x28e9fa9eu,
};

const peace_sm2_JACOBIAN_POINT _peace_sm2_G = {
	{
		0x334c74c7u,
		0x715a4589u,
		0xf2660be1u,
		0x8fe30bbfu,
		0x6a39c994u,
		0x5f990446u,
		0x1f198119u,
		0x32c4ae2cu,
	},
	{
		0x2139f0a0u,
		0x02df32e5u,
		0xc62a4740u,
		0xd0a9877cu,
		0x6b692153u,
		0x59bdcee3u,
		0xf4f6779cu,
		0xbc3736a2u,
	},
	{
		1u,
		0u,
		0u,
		0u,
		0u,
		0u,
		0u,
		0u,
	},
};
const peace_sm2_JACOBIAN_POINT *peace_sm2_G = &_peace_sm2_G;

const peace_sm2_BN peace_sm2_N = {
	0x39d54123u,
	0x53bbf409u,
	0x21c6052bu,
	0x7203df6bu,
	0xffffffffu,
	0xffffffffu,
	0xffffffffu,
	0xfffffffeu,
};

const peace_sm2_BN peace_sm2_U_PLUS_ONE = {
	0x00000000u,
	0x40000000u,
	0xc0000000u,
	0xffffffffu,
	0xffffffffu,
	0xffffffffu,
	0xbfffffffu,
	0x3fffffffu,
};

const peace_sm2_BN peace_sm2_ONE = {1u, 0u, 0u, 0u, 0u, 0u, 0u, 0u};
const peace_sm2_BN peace_sm2_TWO = {2u, 0u, 0u, 0u, 0u, 0u, 0u, 0u};
const peace_sm2_BN peace_sm2_THREE = {3u, 0u, 0u, 0u, 0u, 0u, 0u, 0u};

int peace_sm2_bn_check(const peace_sm2_BN a)
{
	int err = 0;
	int i;
	for (i = 0; i < 8; i++)
	{
		if (a[i] > 0xffffffffu)
		{
			fprintf(stderr, "%s %d: error\n", __FILE__, __LINE__);
			err++;
		}
	}
	if (err)
		return E_ERROR;
	else
		return 1;
}

int peace_sm2_bn_is_zero(const peace_sm2_BN a)
{
	int i;
	for (i = 0; i < 8; i++)
	{
		if (a[i] != 0)
			return 0;
	}
	return 1;
}

int peace_sm2_bn_is_one(const peace_sm2_BN a)
{
	int i;
	if (a[0] != 1)
		return 0;
	for (i = 1; i < 8; i++)
	{
		if (a[i] != 0)
			return 0;
	}
	return 1;
}

/**
 * @brief Converts a big number to a byte array.
 *
 * This function converts the given big number `a` into a 32-byte array `out`.
 * The conversion is done in big-endian order.
 *
 * @param[in] a   Big number to be converted.
 * @param[out] out Output byte array of length 32.
 */
void peace_sm2_bn_to_bytes(const peace_sm2_BN a, uint8_t out[32])
{
    int i;
    for (i = 7; i >= 0; i--)
    {
        uint32_t ai = (uint32_t)a[i];
        PUTU32(out, ai);
        out += sizeof(uint32_t);
    }
}

/**
 * @brief Converts a byte array to a big number.
 *
 * This function converts the given 32-byte array `in` into a big number `r`.
 * The conversion is done in big-endian order.
 *
 * @param[out] r  Big number to store the result.
 * @param[in] in  Input byte array of length 32.
 */
void peace_sm2_bn_from_bytes(peace_sm2_BN r, const uint8_t in[32])
{
    int i;
    for (i = 7; i >= 0; i--)
    {
        r[i] = GETU32(in);
        in += sizeof(uint32_t);
    }
}

/**
 * @brief Converts a hexadecimal character to an integer value.
 *
 * This function converts a single hexadecimal character to its corresponding integer value.
 *
 * @param[in] c   Hexadecimal character ('0'-'9', 'a'-'f', 'A'-'F').
 * @return Integer value of the hexadecimal character or E_ERROR on invalid input.
 */
static int hexchar2int(char c)
{
    if ('0' <= c && c <= '9')
        return c - '0';
    else if ('a' <= c && c <= 'f')
        return c - 'a' + 10;
    else if ('A' <= c && c <= 'F')
        return c - 'A' + 10;
    else
        return E_ERROR;
}

/**
 * @brief Converts a hexadecimal string to a binary byte array.
 *
 * This function converts a hexadecimal string of even length to a binary byte array.
 *
 * @param[in] in  Input hexadecimal string.
 * @param[in] inlen Length of the input string.
 * @param[out] out Output byte array.
 * @return 1 on success, E_ERROR on failure.
 */
static int hex2bin(const char *in, size_t inlen, uint8_t *out)
{
    int c;
    if (inlen % 2)
        return E_ERROR;

    while (inlen)
    {
        if ((c = hexchar2int(*in++)) < 0)
            return E_ERROR;
        *out = (uint8_t)c << 4;
        if ((c = hexchar2int(*in++)) < 0)
            return E_ERROR;
        *out |= (uint8_t)c;
        inlen -= 2;
        out++;
    }
    return 1;
}

/**
 * @brief Converts a big number to a hexadecimal string.
 *
 * This function converts the given big number `a` into a hexadecimal string `hex`.
 *
 * @param[in] a   Big number to be converted.
 * @param[out] hex Output hexadecimal string of length 64.
 */
void peace_sm2_bn_to_hex(const peace_sm2_BN a, char hex[64])
{
    int i;
    for (i = 7; i >= 0; i--)
    {
        int len;
        len = sprintf(hex, "%08x", (uint32_t)a[i]);
        assert(len == 8);
        hex += 8;
    }
}

/**
 * @brief Converts a hexadecimal string to a big number.
 *
 * This function converts the given hexadecimal string `hex` into a big number `r`.
 *
 * @param[out] r  Big number to store the result.
 * @param[in] hex Input hexadecimal string of length 64.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_bn_from_hex(peace_sm2_BN r, const char hex[64])
{
    uint8_t buf[32];
    if (hex2bin(hex, 64, buf) < 0)
        return E_ERROR;
    peace_sm2_bn_from_bytes(r, buf);
    return 1;
}

/**
 * @brief Converts an ASN.1 INTEGER to a big number.
 *
 * This function converts the given ASN.1 INTEGER data `d` into a big number `r`.
 *
 * @param[out] r  Big number to store the result.
 * @param[in] d   Input ASN.1 INTEGER data.
 * @param[in] dlen Length of the input data.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_bn_from_peace_asn1_integer(peace_sm2_BN r, const uint8_t *d, size_t dlen)
{
    uint8_t buf[32] = {0};
    if (!d || dlen == 0)
    {
        return E_ERROR;
    }
    if (dlen > sizeof(buf))
    {
        return E_ERROR;
    }
    memcpy(buf + sizeof(buf) - dlen, d, dlen);
    peace_sm2_bn_from_bytes(r, buf);
    return 1;
}

/**
 * @brief Converts a big number to a binary string of bits.
 *
 * This function converts the given big number `a` into a binary string `bits`.
 *
 * @param[in] a   Big number to be converted.
 * @param[out] bits Output binary string of length 256.
 */
void peace_sm2_bn_to_bits(const peace_sm2_BN a, char bits[256])
{
    int i, j;
    uint64_t w;
    for (i = 7; i >= 0; i--)
    {
        w = a[i];
        for (j = 0; j < 32; j++)
        {
            *bits++ = (w & 0x80000000) ? '1' : '0';
            w <<= 1;
        }
    }
}

/**
 * @brief Compares two big numbers.
 *
 * This function compares two big numbers `a` and `b`.
 *
 * @param[in] a   First big number.
 * @param[in] b   Second big number.
 * @return 1 if `a` is greater than `b`, E_ERROR if `a` is less than `b`, 0 if they are equal.
 */
int peace_sm2_bn_cmp(const peace_sm2_BN a, const peace_sm2_BN b)
{
    int i;
    for (i = 7; i >= 0; i--)
    {
        if (a[i] > b[i])
            return 1;
        if (a[i] < b[i])
            return E_ERROR;
    }
    return 0;
}

/**
 * @brief Checks if a big number is odd.
 *
 * This function checks if the given big number `a` is odd.
 *
 * @param[in] a   Big number to check.
 * @return Non-zero if `a` is odd, zero otherwise.
 */
int peace_sm2_bn_is_odd(const peace_sm2_BN a)
{
    return a[0] & 0x01;
}

/**
 * @brief Sets a big number to a word value.
 *
 * This function sets the given big number `r` to the specified word value `a`.
 *
 * @param[out] r  Big number to set.
 * @param[in] a   Word value to set.
 */
void peace_sm2_bn_set_word(peace_sm2_BN r, uint32_t a)
{
    int i;
    r[0] = a;
    for (i = 1; i < 8; i++)
    {
        r[i] = 0;
    }
}

/**
 * @brief Performs a right shift operation on a big number.
 *
 * This function performs a right shift operation on the given big number `a` by `nbits` bits.
 *
 * @param[out] ret Resulting big number after the shift.
 * @param[in] a   Big number to shift.
 * @param[in] nbits Number of bits to shift.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_bn_rshift(peace_sm2_BN ret, const peace_sm2_BN a, unsigned int nbits)
{
    peace_sm2_BN r;
    int i;

    if (nbits > 31)
    {
        return E_ERROR;
    }
    if (nbits == 0)
    {
        memcpy(ret, a, peace_sm2_sizeof_bn);
    }

    for (i = 0; i < 7; i++)
    {
        r[i] = a[i] >> nbits;
        r[i] |= (a[i + 1] << (32 - nbits)) & 0xffffffffu;
    }
    r[i] = a[i] >> nbits;
    memcpy(ret, r, peace_sm2_sizeof_bn);
    return 1;
}

/**
 * @brief Adds two big numbers.
 *
 * This function adds two big numbers `a` and `b` and stores the result in `r`.
 *
 * @param[out] r  Resulting big number.
 * @param[in] a   First big number.
 * @param[in] b   Second big number.
 */
void peace_sm2_bn_add(peace_sm2_BN r, const peace_sm2_BN a, const peace_sm2_BN b)
{
    int i;
    r[0] = a[0] + b[0];

    for (i = 1; i < 8; i++)
    {
        r[i] = a[i] + b[i] + (r[i - 1] >> 32);
    }
    for (i = 0; i < 7; i++)
    {
        r[i] &= 0xffffffffu;
    }
}

/**
 * @brief Subtracts two big numbers.
 *
 * This function subtracts big number `b` from `a` and stores the result in `ret`.
 *
 * @param[out] ret Resulting big number.
 * @param[in] a   First big number.
 * @param[in] b   Second big number.
 */
void peace_sm2_bn_sub(peace_sm2_BN ret, const peace_sm2_BN a, const peace_sm2_BN b)
{
    int i;
    peace_sm2_BN r;
    r[0] = ((uint64_t)1 << 32) + a[0] - b[0];
    for (i = 1; i < 7; i++)
    {
        r[i] = 0xffffffffu + a[i] - b[i] + (r[i - 1] >> 32);
        r[i - 1] &= 0xffffffffu;
    }
    r[i] = a[i] - b[i] + (r[i - 1] >> 32) - 1;
    r[i - 1] &= 0xffffffffu;
    memcpy(ret, r, peace_sm2_sizeof_bn);
}

/**
 * @brief Generates a random big number within a specified range.
 *
 * This function generates a random big number `r` that is less than `range`.
 *
 * @param[out] r  Resulting random big number.
 * @param[in] range Upper bound for the random number.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_bn_rand_range(peace_sm2_BN r, const peace_sm2_BN range)
{
    uint8_t buf[32];
    do
    {
        srand(time((void *)0));
        for (int a = 0; a < sizeof(buf); a++)
        {
            buf[a] = rand() % 255;
        }

        peace_sm2_bn_from_bytes(r, buf);
    } while (peace_sm2_bn_cmp(r, range) >= 0);
    return 1;
}

/**
 * @brief Adds two finite field elements.
 *
 * This function adds two finite field elements `a` and `b` and stores the result in `r`.
 *
 * @param[out] r  Resulting finite field element.
 * @param[in] a   First finite field element.
 * @param[in] b   Second finite field element.
 */
void peace_sm2_fp_add(peace_sm2_Fp r, const peace_sm2_Fp a, const peace_sm2_Fp b)
{
    peace_sm2_bn_add(r, a, b);
    if (peace_sm2_bn_cmp(r, peace_sm2_P) >= 0)
    {
        peace_sm2_bn_sub(r, r, peace_sm2_P);
    }
}

/**
 * @brief Subtracts two finite field elements.
 *
 * This function subtracts finite field element `b` from `a` and stores the result in `r`.
 *
 * @param[out] r  Resulting finite field element.
 * @param[in] a   First finite field element.
 * @param[in] b   Second finite field element.
 */
void peace_sm2_fp_sub(peace_sm2_Fp r, const peace_sm2_Fp a, const peace_sm2_Fp b)
{
    if (peace_sm2_bn_cmp(a, b) >= 0)
    {
        peace_sm2_bn_sub(r, a, b);
    }
    else
    {
        peace_sm2_BN t;
        peace_sm2_bn_sub(t, peace_sm2_P, b);
        peace_sm2_bn_add(r, t, a);
    }
}

/**
 * @brief Doubles a finite field element.
 *
 * This function doubles the given finite field element `a` and stores the result in `r`.
 *
 * @param[out] r  Resulting finite field element.
 * @param[in] a   Finite field element to double.
 */
void peace_sm2_fp_dbl(peace_sm2_Fp r, const peace_sm2_Fp a)
{
    peace_sm2_fp_add(r, a, a);
}

/**
 * @brief Triples a finite field element.
 *
 * This function triples the given finite field element `a` and stores the result in `r`.
 *
 * @param[out] r  Resulting finite field element.
 * @param[in] a   Finite field element to triple.
 */
void peace_sm2_fp_tri(peace_sm2_Fp r, const peace_sm2_Fp a)
{
    peace_sm2_BN t;
    peace_sm2_fp_dbl(t, a);
    peace_sm2_fp_add(r, t, a);
}

/**
 * @brief Divides a finite field element by 2.
 *
 * This function divides the given finite field element `a` by 2 and stores the result in `r`.
 *
 * @param[out] r  Resulting finite field element.
 * @param[in] a   Finite field element to divide.
 */
void peace_sm2_fp_div2(peace_sm2_Fp r, const peace_sm2_Fp a)
{
    int i;
    memcpy(r, a, peace_sm2_sizeof_bn);
    if (r[0] & 0x01)
    {
        peace_sm2_bn_add(r, r, peace_sm2_P);
    }
    for (i = 0; i < 7; i++)
    {
        r[i] = (r[i] >> 1) | ((r[i + 1] & 0x01) << 31);
    }
    r[i] >>= 1;
}

/**
 * @brief Negates a finite field element.
 *
 * This function negates the given finite field element `a` and stores the result in `r`.
 *
 * @param[out] r  Resulting finite field element.
 * @param[in] a   Finite field element to negate.
 */
void peace_sm2_fp_neg(peace_sm2_Fp r, const peace_sm2_Fp a)
{
    if (peace_sm2_bn_is_zero(a))
    {
        memcpy(r, a, peace_sm2_sizeof_bn);
    }
    else
    {
        peace_sm2_bn_sub(r, peace_sm2_P, a);
    }
}

/**
 * @brief Multiplies two finite field elements.
 *
 * This function multiplies two finite field elements `a` and `b` and stores the result in `r`.
 *
 * @param[out] r  Resulting finite field element.
 * @param[in] a   First finite field element.
 * @param[in] b   Second finite field element.
 */

void peace_sm2_fp_mul(peace_sm2_Fp r, const peace_sm2_Fp a, const peace_sm2_Fp b)
{
	int i, j;
	uint64_t s[16] = {0};
	peace_sm2_BN d = {0};
	uint64_t u;


	for (i = 0; i < 8; i++)
	{
		u = 0;
		for (j = 0; j < 8; j++)
		{
			u = s[i + j] + a[i] * b[j] + u;
			s[i + j] = u & 0xffffffffu;
			u >>= 32;
		}
		s[i + 8] = u;
	}

	r[0] = s[0] + s[8] + s[9] + s[10] + s[11] + s[12] + ((s[13] + s[14] + s[15]) << 1);
	r[1] = s[1] + s[9] + s[10] + s[11] + s[12] + s[13] + ((s[14] + s[15]) << 1);
	r[2] = s[2];
	r[3] = s[3] + s[8] + s[11] + s[12] + s[14] + s[15] + (s[13] << 1);
	r[4] = s[4] + s[9] + s[12] + s[13] + s[15] + (s[14] << 1);
	r[5] = s[5] + s[10] + s[13] + s[14] + (s[15] << 1);
	r[6] = s[6] + s[11] + s[14] + s[15];
	r[7] = s[7] + s[8] + s[9] + s[10] + s[11] + s[15] + ((s[12] + s[13] + s[14] + s[15]) << 1);

	for (i = 1; i < 8; i++)
	{
		r[i] += r[i - 1] >> 32;
		r[i - 1] &= 0xffffffffu;
	}

	d[2] = s[8] + s[9] + s[13] + s[14];
	d[3] = d[2] >> 32;
	d[2] &= 0xffffffffu;
	peace_sm2_bn_sub(r, r, d);

	while (peace_sm2_bn_cmp(r, peace_sm2_P) >= 0)
	{
		peace_sm2_bn_sub(r, r, peace_sm2_P);
	}
}

/**
 * @brief Squares the input and stores the result in r.
 * @param[in,out] r The result of squaring a, passed as an output parameter.
 * @param[in] a The input to be squared.
 */
void peace_sm2_fp_sqr(peace_sm2_Fp r, const peace_sm2_Fp a)
{
    // Multiply a by itself and store the result in r
    peace_sm2_fp_mul(r, a, a);
}

/**
 * @brief Computes the exponentiation of a to the power e and stores the result in r.
 * @param[in,out] r The result of a^e, passed as an output parameter.
 * @param[in] a The base of the exponentiation.
 * @param[in] e The exponent.
 */
void peace_sm2_fp_exp(peace_sm2_Fp r, const peace_sm2_Fp a, const peace_sm2_Fp e)
{
    peace_sm2_BN t;
    uint32_t w;
    int i, j;

    // Initialize t to 1
    peace_sm2_bn_set_word(t, 1);
    for (i = 7; i >= 0; i--)
    {
        w = (uint32_t)e[i];
        for (j = 0; j < 32; j++)
        {
            // Square t
            peace_sm2_fp_sqr(t, t);
            if (w & 0x80000000)
                // If the current bit of e is 1, multiply t by a
                peace_sm2_fp_mul(t, t, a);
            w <<= 1;
        }
    }

    // Copy the result from t to r
    memcpy(r, t, peace_sm2_sizeof_bn);
}

/**
 * @brief Computes the modular inverse of a and stores the result in r.
 * @param[in,out] r The result of the modular inverse of a, passed as an output parameter.
 * @param[in] a The input for which the modular inverse is computed.
 */
void peace_sm2_fp_inv(peace_sm2_Fp r, const peace_sm2_Fp a)
{
    peace_sm2_BN a1, a2, a3, a4, a5;
    int i;

    // Intermediate calculations for the modular inverse
    peace_sm2_fp_sqr(a1, a);
    peace_sm2_fp_mul(a2, a1, a);
    peace_sm2_fp_sqr(a3, a2);
    peace_sm2_fp_sqr(a3, a3);
    peace_sm2_fp_mul(a3, a3, a2);
    peace_sm2_fp_sqr(a4, a3);
    peace_sm2_fp_sqr(a4, a4);
    peace_sm2_fp_sqr(a4, a4);
    peace_sm2_fp_sqr(a4, a4);
    peace_sm2_fp_mul(a4, a4, a3);
    peace_sm2_fp_sqr(a5, a4);
    for (i = 1; i < 8; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a5, a5, a4);
    for (i = 0; i < 8; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a5, a5, a4);
    for (i = 0; i < 4; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a5, a5, a3);
    peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a5, a5, a2);
    peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a5, a5, a);
    peace_sm2_fp_sqr(a4, a5);
    peace_sm2_fp_mul(a3, a4, a1);
    peace_sm2_fp_sqr(a5, a4);
    for (i = 1; i < 31; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a4, a5, a4);
    peace_sm2_fp_sqr(a4, a4);
    peace_sm2_fp_mul(a4, a4, a);
    peace_sm2_fp_mul(a3, a4, a2);
    for (i = 0; i < 33; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a2, a5, a3);
    peace_sm2_fp_mul(a3, a2, a3);
    for (i = 0; i < 32; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a2, a5, a3);
    peace_sm2_fp_mul(a3, a2, a3);
    peace_sm2_fp_mul(a4, a2, a4);
    for (i = 0; i < 32; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a2, a5, a3);
    peace_sm2_fp_mul(a3, a2, a3);
    peace_sm2_fp_mul(a4, a2, a4);
    for (i = 0; i < 32; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a2, a5, a3);
    peace_sm2_fp_mul(a3, a2, a3);
    peace_sm2_fp_mul(a4, a2, a4);
    for (i = 0; i < 32; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(a2, a5, a3);
    peace_sm2_fp_mul(a3, a2, a3);
    peace_sm2_fp_mul(a4, a2, a4);
    for (i = 0; i < 32; i++)
        peace_sm2_fp_sqr(a5, a5);
    peace_sm2_fp_mul(r, a4, a5);

    // Clear intermediate variables
    memset(a1, 0, peace_sm2_sizeof_bn);
    memset(a2, 0, peace_sm2_sizeof_bn);
    memset(a3, 0, peace_sm2_sizeof_bn);
    memset(a4, 0, peace_sm2_sizeof_bn);
    memset(a5, 0, peace_sm2_sizeof_bn);
}

/**
 * @brief Adds two numbers modulo N and stores the result in r.
 * @param[in,out] r The result of (a + b) % N, passed as an output parameter.
 * @param[in] a The first operand.
 * @param[in] b The second operand.
 */
void peace_sm2_fn_add(peace_sm2_Fn r, const peace_sm2_Fn a, const peace_sm2_Fn b)
{
    // Add a and b
    peace_sm2_bn_add(r, a, b);
    if (peace_sm2_bn_cmp(r, peace_sm2_N) >= 0)
    {
        // If the result is greater than or equal to N, subtract N
        peace_sm2_bn_sub(r, r, peace_sm2_N);
    }
}

/**
 * @brief Subtracts b from a modulo N and stores the result in r.
 * @param[in,out] r The result of (a - b) % N, passed as an output parameter.
 * @param[in] a The minuend.
 * @param[in] b The subtrahend.
 */
void peace_sm2_fn_sub(peace_sm2_Fn r, const peace_sm2_Fn a, const peace_sm2_Fn b)
{
    if (peace_sm2_bn_cmp(a, b) >= 0)
    {
        // If a is greater than or equal to b, subtract b from a
        peace_sm2_bn_sub(r, a, b);
    }
    else
    {
        peace_sm2_BN t;
        // If a is less than b, add N to a and then subtract b
        peace_sm2_bn_add(t, a, peace_sm2_N);
        peace_sm2_bn_sub(r, t, b);
    }
}

/**
 * @brief Computes the negative of a modulo N and stores the result in r.
 * @param[in,out] r The result of (-a) % N, passed as an output parameter.
 * @param[in] a The input for which the negative is computed.
 */
void peace_sm2_fn_neg(peace_sm2_Fn r, const peace_sm2_Fn a)
{
    if (peace_sm2_bn_is_zero(a))
    {
        // If a is zero, copy a to r
        memcpy(r, a, peace_sm2_sizeof_bn);
    }
    else
    {
        // Compute N - a
        peace_sm2_bn_sub(r, peace_sm2_N, a);
    }
}

/**
 * @brief Compares two 288-bit numbers.
 * @param[in] a The first 288-bit number.
 * @param[in] b The second 288-bit number.
 * @return 1 if a > b, -1 if a < b, 0 if a == b.
 */
static int peace_sm2_bn288_cmp(const uint64_t a[9], const uint64_t b[9])
{
    int i;
    for (i = 8; i >= 0; i--)
    {
        if (a[i] > b[i])
            return 1;
        if (a[i] < b[i])
            return -1;
    }
    return 0;
}

/**
 * @brief Adds two 288-bit numbers and stores the result in r.
 * @param[in,out] r The result of a + b, passed as an output parameter.
 * @param[in] a The first 288-bit number.
 * @param[in] b The second 288-bit number.
 */
static void peace_sm2_bn288_add(uint64_t r[9], const uint64_t a[9], const uint64_t b[9])
{
    int i;
    r[0] = a[0] + b[0];
    for (i = 1; i < 9; i++)
    {
        r[i] = a[i] + b[i] + (r[i - 1] >> 32);
    }
    for (i = 0; i < 8; i++)
    {
        r[i] &= 0xffffffffu;
    }
}

/**
 * @brief Subtracts two 288-bit numbers and stores the result in ret.
 * @param[in,out] ret The result of a - b, passed as an output parameter.
 * @param[in] a The minuend.
 * @param[in] b The subtrahend.
 */
static void peace_sm2_bn288_sub(uint64_t ret[9], const uint64_t a[9], const uint64_t b[9])
{
    int i;
    uint64_t r[9];

    r[0] = ((uint64_t)1 << 32) + a[0] - b[0];
    for (i = 1; i < 8; i++)
    {
        r[i] = 0xffffffffu + a[i] - b[i] + (r[i - 1] >> 32);
        r[i - 1] &= 0xffffffffu;
    }
    r[i] = a[i] - b[i] + (r[i - 1] >> 32) - 1;
    r[i - 1] &= 0xffffffffu;

    for (i = 0; i < 9; i++)
    {
        ret[i] = r[i];
    }
}

/**
 * @brief Multiplies two big numbers and returns the result.
 * 
 * This function performs modular multiplication of two large numbers (a, b) and stores the result in ret.
 * The operation is performed modulo n, where n is a predefined constant.
 * 
 * @param [out] ret Pointer to the resulting big number.
 * @param [in] a Pointer to the first operand.
 * @param [in] b Pointer to the second operand.
 */
void peace_sm2_fn_mul(peace_sm2_BN ret, const peace_sm2_BN a, const peace_sm2_BN b)
{
    // Temporary variables for intermediate results
    peace_sm2_BN r;
    static const uint64_t mu[9] = {
        0xf15149a0,
        0x12ac6361,
        0xfa323c01,
        0x8dfc2096,
        1,
        1,
        1,
        1,
        1,
    };

    uint64_t s[18];
    uint64_t zh[9];
    uint64_t zl[9];
    uint64_t q[9];
    uint64_t w;
    int i, j;

    /* Calculate the product z = a * b */
    for (i = 0; i < 8; i++)
    {
        s[i] = 0;
    }
    for (i = 0; i < 8; i++)
    {
        w = 0;
        for (j = 0; j < 8; j++)
        {
            w += s[i + j] + a[i] * b[j];
            s[i + j] = w & 0xffffffffu;
            w >>= 32;
        }
        s[i + 8] = w;
    }

    // Split the product into high and low parts
    for (i = 0; i < 9; i++)
    {
        zl[i] = s[i];
        zh[i] = s[7 + i];
    }

    // Calculate q = zh * mu
    for (i = 0; i < 9; i++)
    {
        s[i] = 0;
    }
    for (i = 0; i < 9; i++)
    {
        w = 0;
        for (j = 0; j < 9; j++)
        {
            w += s[i + j] + zh[i] * mu[j];
            s[i + j] = w & 0xffffffffu;
            w >>= 32;
        }
        s[i + 9] = w;
    }
    for (i = 0; i < 8; i++)
    {
        q[i] = s[9 + i];
    }

    // Reduce q by n
    for (i = 0; i < 17; i++)
    {
        s[i] = 0;
    }
    for (i = 0; i < 8; i++)
    {
        w = 0;
        for (j = 0; j < 8; j++)
        {
            w += s[i + j] + q[i] * peace_sm2_N[j];
            s[i + j] = w & 0xffffffffu;
            w >>= 32;
        }
        s[i + 8] = w;
    }
    for (i = 0; i < 9; i++)
    {
        q[i] = s[i];
    }

    // Calculate r = zl - q (mod (2^32)^9)
    if (peace_sm2_bn288_cmp(zl, q))
    {
        peace_sm2_bn288_sub(zl, zl, q);
    }
    else
    {
        uint64_t c[9] = {0, 0, 0, 0, 0, 0, 0, 0, 0x100000000};
        peace_sm2_bn288_sub(q, c, q);
        peace_sm2_bn288_add(zl, q, zl);
    }

    for (i = 0; i < 8; i++)
    {
        r[i] = zl[i];
    }
    r[7] += zl[8] << 32;

    // Ensure r is less than n
    while (peace_sm2_bn_cmp(r, peace_sm2_N) >= 0)
    {
        peace_sm2_bn_sub(r, r, peace_sm2_N);
    }

    // Copy the result to the output
    memcpy(ret, r, peace_sm2_sizeof_bn);
}

/**
 * @brief Multiplies a big number by a word and returns the result.
 * 
 * This function multiplies a big number (a) by a 32-bit word (b) and stores the result in r.
 * 
 * @param [out] r Pointer to the resulting big number.
 * @param [in] a Pointer to the operand.
 * @param [in] b 32-bit word to multiply with.
 */
void peace_sm2_fn_mul_word(peace_sm2_Fn r, const peace_sm2_Fn a, uint32_t b)
{
    peace_sm2_Fn t;
    peace_sm2_bn_set_word(t, b);
    peace_sm2_fn_mul(r, a, t);
}

/**
 * @brief Squares a big number and returns the result.
 * 
 * This function squares a big number (a) and stores the result in r.
 * 
 * @param [out] r Pointer to the resulting big number.
 * @param [in] a Pointer to the operand.
 */
void peace_sm2_fn_sqr(peace_sm2_BN r, const peace_sm2_BN a)
{
    peace_sm2_fn_mul(r, a, a);
}

/**
 * @brief Raises a big number to a power and returns the result.
 * 
 * This function raises a big number (a) to the power of another big number (e) and stores the result in r.
 * 
 * @param [out] r Pointer to the resulting big number.
 * @param [in] a Pointer to the base.
 * @param [in] e Pointer to the exponent.
 */
void peace_sm2_fn_exp(peace_sm2_BN r, const peace_sm2_BN a, const peace_sm2_BN e)
{
    peace_sm2_BN t;
    uint32_t w;
    int i, j;

    peace_sm2_bn_set_word(t, 1);
    for (i = 7; i >= 0; i--)
    {
        w = (uint32_t)e[i];
        for (j = 0; j < 32; j++)
        {
            peace_sm2_fn_sqr(t, t);
            if (w & 0x80000000)
            {
                peace_sm2_fn_mul(t, t, a);
            }
            w <<= 1;
        }
    }
    memcpy(r, t, peace_sm2_sizeof_bn);
}

/**
 * @brief Calculates the modular inverse of a big number and returns the result.
 * 
 * This function calculates the modular inverse of a big number (a) and stores the result in r.
 * 
 * @param [out] r Pointer to the resulting big number.
 * @param [in] a Pointer to the operand.
 */
void peace_sm2_fn_inv(peace_sm2_BN r, const peace_sm2_BN a)
{
    peace_sm2_BN e;
    peace_sm2_bn_sub(e, peace_sm2_N, peace_sm2_TWO);
    peace_sm2_fn_exp(r, a, e);
}

/**
 * @brief Generates a random big number within a specified range.
 * 
 * This function generates a random big number (r) within the range of 0 to n-1.
 * 
 * @param [out] r Pointer to the resulting big number.
 * @return 1 if successful, E_ERROR otherwise.
 */
int peace_sm2_fn_rand(peace_sm2_BN r)
{
    if (peace_sm2_bn_rand_range(r, peace_sm2_N) != 1)
    {
        return E_ERROR;
    }
    return 1;
}

/**
 * @brief Initializes a Jacobian point structure.
 * 
 * This function initializes a Jacobian point structure (R) to a default value.
 * 
 * @param [out] R Pointer to the Jacobian point structure.
 */
void peace_sm2_jacobian_point_init(peace_sm2_JACOBIAN_POINT *R)
{
    memset(R, 0, sizeof(peace_sm2_JACOBIAN_POINT));
    R->X[0] = 1;
    R->Y[0] = 1;
}

/**
 * @brief Checks if a Jacobian point is at infinity.
 * 
 * This function checks if a Jacobian point (P) is at infinity.
 * 
 * @param [in] P Pointer to the Jacobian point.
 * @return 1 if the point is at infinity, 0 otherwise.
 */
int peace_sm2_jacobian_point_is_at_infinity(const peace_sm2_JACOBIAN_POINT *P)
{
    return peace_sm2_bn_is_zero(P->Z);
}

/**
 * @brief Sets the x and y coordinates of a Jacobian point.
 * 
 * This function sets the x and y coordinates of a Jacobian point (R) to the given values.
 * 
 * @param [out] R Pointer to the Jacobian point.
 * @param [in] x Pointer to the x-coordinate.
 * @param [in] y Pointer to the y-coordinate.
 */
void peace_sm2_jacobian_point_set_xy(peace_sm2_JACOBIAN_POINT *R, const peace_sm2_BN x, const peace_sm2_BN y)
{
    memcpy(R->X, x, peace_sm2_sizeof_bn);
    memcpy(R->Y, y, peace_sm2_sizeof_bn);
    peace_sm2_bn_set_word(R->Z, 1);
}


/**
 * @brief Converts a Jacobian point to affine coordinates.
 *
 * This function converts a point from Jacobian coordinates to affine coordinates.
 * If the point is at infinity, the function will handle it accordingly.
 *
 * @param P Pointer to the input Jacobian point.
 * @param x Pointer to the output affine x-coordinate.
 * @param y Pointer to the output affine y-coordinate (optional).
 */
void peace_sm2_jacobian_point_get_xy(const peace_sm2_JACOBIAN_POINT *P, peace_sm2_BN x, peace_sm2_BN y)
{
    if (peace_sm2_bn_is_one(P->Z))
    {
        memcpy(x, P->X, peace_sm2_sizeof_bn);
        if (y)
        {
            memcpy(y, P->Y, peace_sm2_sizeof_bn);
        }
    }
    else
    {
        peace_sm2_BN z_inv;
        peace_sm2_fp_inv(z_inv, P->Z);
        if (y)
        {
            peace_sm2_fp_mul(y, P->Y, z_inv);
        }
        peace_sm2_fp_sqr(z_inv, z_inv);
        peace_sm2_fp_mul(x, P->X, z_inv);
        if (y)
        {
            peace_sm2_fp_mul(y, y, z_inv);
        }
    }
}

/**
 * @brief Checks if a Jacobian point lies on the SM2 elliptic curve.
 *
 * This function verifies whether the given Jacobian point lies on the SM2 elliptic curve.
 * It performs the necessary calculations to check the curve equation.
 *
 * @param P Pointer to the Jacobian point to be checked.
 * @return 1 if the point is on the curve, E_ERROR otherwise.
 */
int peace_sm2_jacobian_point_is_on_curve(const peace_sm2_JACOBIAN_POINT *P)
{
    peace_sm2_BN t0;
    peace_sm2_BN t1;
    peace_sm2_BN t2;

    if (peace_sm2_bn_is_one(P->Z))
    {
        peace_sm2_fp_sqr(t0, P->Y);
        peace_sm2_fp_add(t0, t0, P->X);
        peace_sm2_fp_add(t0, t0, P->X);
        peace_sm2_fp_add(t0, t0, P->X);
        peace_sm2_fp_sqr(t1, P->X);
        peace_sm2_fp_mul(t1, t1, P->X);
        peace_sm2_fp_add(t1, t1, peace_sm2_B);
    }
    else
    {
        peace_sm2_fp_sqr(t0, P->Y);
        peace_sm2_fp_sqr(t1, P->Z);
        peace_sm2_fp_sqr(t2, t1);
        peace_sm2_fp_mul(t1, t1, t2);
        peace_sm2_fp_mul(t1, t1, peace_sm2_B);
        peace_sm2_fp_mul(t2, t2, P->X);
        peace_sm2_fp_add(t0, t0, t2);
        peace_sm2_fp_add(t0, t0, t2);
        peace_sm2_fp_add(t0, t0, t2);
        peace_sm2_fp_sqr(t2, P->X);
        peace_sm2_fp_mul(t2, t2, P->X);
        peace_sm2_fp_add(t1, t1, t2);
    }

    if (peace_sm2_bn_cmp(t0, t1) != 0)
    {
        return E_ERROR;
    }
    return 1;
}

/**
 * @brief Negates a Jacobian point.
 *
 * This function negates the y-coordinate of the given Jacobian point.
 * The x and z coordinates remain unchanged.
 *
 * @param R Pointer to the output Jacobian point.
 * @param P Pointer to the input Jacobian point.
 */
void peace_sm2_jacobian_point_neg(peace_sm2_JACOBIAN_POINT *R, const peace_sm2_JACOBIAN_POINT *P)
{
    memcpy(R->X, P->X, peace_sm2_sizeof_bn);
    peace_sm2_fp_neg(R->Y, P->Y);
    memcpy(R->Z, P->Z, peace_sm2_sizeof_bn);
}

/**
 * @brief Doubles a Jacobian point.
 *
 * This function doubles the given Jacobian point using the addition formula for elliptic curves.
 * If the point is at infinity, it returns the point itself.
 *
 * @param R Pointer to the output Jacobian point.
 * @param P Pointer to the input Jacobian point.
 */
void peace_sm2_jacobian_point_dbl(peace_sm2_JACOBIAN_POINT *R, const peace_sm2_JACOBIAN_POINT *P)
{
    const uint64_t *X1 = P->X;
    const uint64_t *Y1 = P->Y;
    const uint64_t *Z1 = P->Z;
    peace_sm2_BN T1;
    peace_sm2_BN T2;
    peace_sm2_BN T3;
    peace_sm2_BN X3;
    peace_sm2_BN Y3;
    peace_sm2_BN Z3;

    if (peace_sm2_jacobian_point_is_at_infinity(P))
    {
        peace_sm2_jacobian_point_copy(R, P);
        return;
    }

    peace_sm2_fp_sqr(T1, Z1);        
    peace_sm2_fp_sub(T2, X1, T1);
    peace_sm2_fp_add(T1, X1, T1); 
    peace_sm2_fp_mul(T2, T2, T1); 
    peace_sm2_fp_tri(T2, T2);        
    peace_sm2_fp_dbl(Y3, Y1);        
    peace_sm2_fp_mul(Z3, Y3, Z1); 
    peace_sm2_fp_sqr(Y3, Y3);        
    peace_sm2_fp_mul(T3, Y3, X1); 
    peace_sm2_fp_sqr(Y3, Y3);        
    peace_sm2_fp_div2(Y3, Y3);    
    peace_sm2_fp_sqr(X3, T2);        
    peace_sm2_fp_dbl(T1, T3);        
    peace_sm2_fp_sub(X3, X3, T1); 
    peace_sm2_fp_sub(T1, T3, X3); 
    peace_sm2_fp_mul(T1, T1, T2); 
    peace_sm2_fp_sub(Y3, T1, Y3); 

    memcpy(R->X, X3, peace_sm2_sizeof_bn);
    memcpy(R->Y, Y3, peace_sm2_sizeof_bn);
    memcpy(R->Z, Z3, peace_sm2_sizeof_bn);
}

/**
* @brief Adds two Jacobian points.
* @param[out] R Result of the addition.
* @param[in] P First point to add.
* @param[in] Q Second point to add.
* @pre P and Q must be valid Jacobian points.
* @post R will contain the result of the addition.
* @return None
*/
void peace_sm2_jacobian_point_add(peace_sm2_JACOBIAN_POINT *R, const peace_sm2_JACOBIAN_POINT *P, const peace_sm2_JACOBIAN_POINT *Q)
{
    const uint64_t *X1 = P->X;
    const uint64_t *Y1 = P->Y;
    const uint64_t *Z1 = P->Z;
    const uint64_t *x2 = Q->X;
    const uint64_t *y2 = Q->Y;
    peace_sm2_BN T1;
    peace_sm2_BN T2;
    peace_sm2_BN T3;
    peace_sm2_BN T4;
    peace_sm2_BN X3;
    peace_sm2_BN Y3;
    peace_sm2_BN Z3;

    if (peace_sm2_jacobian_point_is_at_infinity(Q))
    {
        peace_sm2_jacobian_point_copy(R, P);
        return;
    }

    if (peace_sm2_jacobian_point_is_at_infinity(P))
    {
        peace_sm2_jacobian_point_copy(R, Q);
        return;
    }

    assert(peace_sm2_bn_is_one(Q->Z));

    peace_sm2_fp_sqr(T1, Z1);
    peace_sm2_fp_mul(T2, T1, Z1);
    peace_sm2_fp_mul(T1, T1, x2);
    peace_sm2_fp_mul(T2, T2, y2);
    peace_sm2_fp_sub(T1, T1, X1);
    peace_sm2_fp_sub(T2, T2, Y1);
    if (peace_sm2_bn_is_zero(T1))
    {
        if (peace_sm2_bn_is_zero(T2))
        {
            peace_sm2_JACOBIAN_POINT _Q, *Q = &_Q;
            peace_sm2_jacobian_point_set_xy(Q, x2, y2);

            peace_sm2_jacobian_point_dbl(R, Q);
            return;
        }
        else
        {
            peace_sm2_jacobian_point_init(R);
            return;
        }
    }
    peace_sm2_fp_mul(Z3, Z1, T1);
    peace_sm2_fp_sqr(T3, T1);
    peace_sm2_fp_mul(T4, T3, T1);
    peace_sm2_fp_mul(T3, T3, X1);
    peace_sm2_fp_dbl(T1, T3);
    peace_sm2_fp_sqr(X3, T2);
    peace_sm2_fp_sub(X3, X3, T1);
    peace_sm2_fp_sub(X3, X3, T4);
    peace_sm2_fp_sub(T3, T3, X3);
    peace_sm2_fp_mul(T3, T3, T2);
    peace_sm2_fp_mul(T4, T4, Y1);
    peace_sm2_fp_sub(Y3, T3, T4);

    memcpy(R->X, X3, peace_sm2_sizeof_bn);
    memcpy(R->Y, Y3, peace_sm2_sizeof_bn);
    memcpy(R->Z, Z3, peace_sm2_sizeof_bn);
}

/**
* @brief Subtracts one Jacobian point from another.
* @param[out] R Result of the subtraction.
* @param[in] P Point to subtract from.
* @param[in] Q Point to subtract.
* @pre P and Q must be valid Jacobian points.
* @post R will contain the result of the subtraction.
* @return None
*/
void peace_sm2_jacobian_point_sub(peace_sm2_JACOBIAN_POINT *R, const peace_sm2_JACOBIAN_POINT *P, const peace_sm2_JACOBIAN_POINT *Q)
{
    peace_sm2_JACOBIAN_POINT _T, *T = &_T;
    peace_sm2_jacobian_point_neg(T, Q);
    peace_sm2_jacobian_point_add(R, P, T);
}

/**
* @brief Multiplies a Jacobian point by a scalar.
* @param[out] R Result of the multiplication.
* @param[in] k Scalar to multiply by.
* @param[in] P Point to multiply.
* @pre P must be a valid Jacobian point.
* @post R will contain the result of the multiplication.
* @return None
*/
void peace_sm2_jacobian_point_mul(peace_sm2_JACOBIAN_POINT *R, const peace_sm2_BN k, const peace_sm2_JACOBIAN_POINT *P)
{
    char bits[257] = {0};
    peace_sm2_JACOBIAN_POINT _Q, *Q = &_Q;
    peace_sm2_JACOBIAN_POINT _T, *T = &_T;
    int i;

    if (!peace_sm2_bn_is_one(P->Z))
    {
        peace_sm2_BN x;
        peace_sm2_BN y;
        peace_sm2_jacobian_point_get_xy(P, x, y);
        peace_sm2_jacobian_point_set_xy(T, x, y);
        P = T;
    }

    peace_sm2_jacobian_point_init(Q);
    peace_sm2_bn_to_bits(k, bits);
    for (i = 0; i < 256; i++)
    {
        peace_sm2_jacobian_point_dbl(Q, Q);
        if (bits[i] == '1')
        {
            peace_sm2_jacobian_point_add(Q, Q, P);
        }
    }
    peace_sm2_jacobian_point_copy(R, Q);
}

/**
* @brief Converts a Jacobian point to bytes.
* @param[in] P Point to convert.
* @param[out] out Buffer to store the result.
* @pre P must be a valid Jacobian point.
* @post out will contain the byte representation of the point.
* @return None
*/
void peace_sm2_jacobian_point_to_bytes(const peace_sm2_JACOBIAN_POINT *P, uint8_t out[64])
{
    peace_sm2_BN x;
    peace_sm2_BN y;
    peace_sm2_jacobian_point_get_xy(P, x, y);
    peace_sm2_bn_to_bytes(x, out);
    peace_sm2_bn_to_bytes(y, out + 32);
}

/**
* @brief Converts bytes to a Jacobian point.
* @param[out] P Point to store the result.
* @param[in] in Buffer containing the byte representation.
* @pre in must be a valid byte representation of a Jacobian point.
* @post P will contain the Jacobian point.
* @return None
*/
void peace_sm2_jacobian_point_from_bytes(peace_sm2_JACOBIAN_POINT *P, const uint8_t in[64])
{
    peace_sm2_bn_from_bytes(P->X, in);
    peace_sm2_bn_from_bytes(P->Y, in + 32);
    peace_sm2_bn_set_word(P->Z, 1);
    /* should we check if peace_sm2_jacobian_point_is_on_curve */
}

/**
* @brief Multiplies the generator point by a scalar.
* @param[out] R Result of the multiplication.
* @param[in] k Scalar to multiply by.
* @pre k must be a valid scalar.
* @post R will contain the result of the multiplication.
* @return None
*/
void peace_sm2_jacobian_point_mul_generator(peace_sm2_JACOBIAN_POINT *R, const peace_sm2_BN k)
{
    peace_sm2_jacobian_point_mul(R, k, peace_sm2_G);
}

/**
* @brief Computes t * P + s * G.
* @param[out] R Result of the computation.
* @param[in] t Scalar to multiply P by.
* @param[in] P Point to multiply.
* @param[in] s Scalar to multiply the generator by.
* @pre P must be a valid Jacobian point.
* @pre t and s must be valid scalars.
* @post R will contain the result of the computation.
* @return None
*/
void peace_sm2_jacobian_point_mul_sum(peace_sm2_JACOBIAN_POINT *R, const peace_sm2_BN t, const peace_sm2_JACOBIAN_POINT *P, const peace_sm2_BN s)
{
    peace_sm2_JACOBIAN_POINT _sG, *sG = &_sG;
    peace_sm2_BN x;
    peace_sm2_BN y;

    peace_sm2_jacobian_point_mul_generator(sG, s);

    peace_sm2_jacobian_point_mul(R, t, P);
    peace_sm2_jacobian_point_get_xy(R, x, y);
    peace_sm2_jacobian_point_set_xy(R, x, y);

    peace_sm2_jacobian_point_add(R, sG, R);
}

/**
* @brief Converts a hexadecimal string to a Jacobian point.
* @param[out] P Point to store the result.
* @param[in] hex Hexadecimal string representation of the point.
* @pre hex must be a valid hexadecimal string.
* @post P will contain the Jacobian point.
* @return None
*/
void peace_sm2_jacobian_point_from_hex(peace_sm2_JACOBIAN_POINT *P, const char hex[64 * 2])
{
    peace_sm2_bn_from_hex(P->X, hex);
    peace_sm2_bn_from_hex(P->Y, hex + 64);
    peace_sm2_bn_set_word(P->Z, 1);
}

/**
* @brief Checks if a Jacobian point is equal to a hexadecimal string.
* @param[in] P Point to compare.
* @param[in] hex Hexadecimal string representation of the point.
* @pre P must be a valid Jacobian point.
* @pre hex must be a valid hexadecimal string.
* @return 1 if the points are equal, 0 otherwise.
*/
int peace_sm2_jacobian_point_equ_hex(const peace_sm2_JACOBIAN_POINT *P, const char hex[128])
{
    peace_sm2_BN x;
    peace_sm2_BN y;
    peace_sm2_JACOBIAN_POINT _T, *T = &_T;

    peace_sm2_jacobian_point_get_xy(P, x, y);
    peace_sm2_jacobian_point_from_hex(T, hex);

    return (peace_sm2_bn_cmp(x, T->X) == 0) && (peace_sm2_bn_cmp(y, T->Y) == 0);
}

/**
* @brief Checks if a point is on the curve.
* @param[in] P Point to check.
* @pre P must be a valid point.
* @return 1 if the point is on the curve, 0 otherwise.
*/
int peace_sm2_point_is_on_curve(const peace_sm2_POINT *P)
{
    peace_sm2_JACOBIAN_POINT T;
    peace_sm2_jacobian_point_from_bytes(&T, (const uint8_t *)P);
    return peace_sm2_jacobian_point_is_on_curve(&T);
}

/**
* @brief Creates a point from its x-coordinate and a parity bit.
* @param[out] P Point to store the result.
* @param[in] x X-coordinate of the point.
* @param[in] y Parity bit of the y-coordinate.
* @pre x must be a valid x-coordinate.
* @post P will contain the point.
* @return 1 if successful, E_ERROR otherwise.
*/
int peace_sm2_point_from_x(peace_sm2_POINT *P, const uint8_t x[32], int y)
{
    peace_sm2_BN _x, _y, _g, _z;
    peace_sm2_bn_from_bytes(_x, x);

    peace_sm2_fp_sqr(_g, _x);
    peace_sm2_fp_sub(_g, _g, peace_sm2_THREE);
    peace_sm2_fp_mul(_g, _g, _x);
    peace_sm2_fp_add(_g, _g, peace_sm2_B);

    peace_sm2_fp_exp(_y, _g, peace_sm2_U_PLUS_ONE);

    peace_sm2_fp_sqr(_z, _y);
    if (peace_sm2_bn_cmp(_z, _g))
    {
        return E_ERROR;
    }

    if ((y == 0x02 && peace_sm2_bn_is_odd(_y)) || ((y == 0x03) && !peace_sm2_bn_is_odd(_y)))
    {
        peace_sm2_fp_neg(_y, _y);
    }

    peace_sm2_bn_to_bytes(_x, P->x);
    peace_sm2_bn_to_bytes(_y, P->y);

    memset(_x, 0, peace_sm2_sizeof_bn);
    memset(_y, 0, peace_sm2_sizeof_bn);
    memset(_g, 0, peace_sm2_sizeof_bn);
    memset(_z, 0, peace_sm2_sizeof_bn);

    if (peace_sm2_point_is_on_curve(P) != 1)
    {
        return E_ERROR;
    }
    return 1;
}

/**
* @brief Creates a point from its x and y coordinates.
* @param[out] P Point to store the result.
* @param[in] x X-coordinate of the point.
* @param[in] y Y-coordinate of the point.
* @pre x and y must be valid coordinates.
* @post P will contain the point.
* @return 1 if successful, E_ERROR otherwise.
*/
int peace_sm2_point_from_xy(peace_sm2_POINT *P, const uint8_t x[32], const uint8_t y[32])
{
    memcpy(P->x, x, 32);
    memcpy(P->y, y, 32);
    return peace_sm2_point_is_on_curve(P);
}

/**
* @brief Adds two points.
* @param[out] R Result of the addition.
* @param[in] P First point to add.
* @param[in] Q Second point to add.
* @pre P and Q must be valid points.
* @post R will contain the result of the addition.
* @return 1 if successful, E_ERROR otherwise.
*/
int peace_sm2_point_add(peace_sm2_POINT *R, const peace_sm2_POINT *P, const peace_sm2_POINT *Q)
{
    peace_sm2_JACOBIAN_POINT P_;
    peace_sm2_JACOBIAN_POINT Q_;

    peace_sm2_jacobian_point_from_bytes(&P_, (uint8_t *)P);
    peace_sm2_jacobian_point_from_bytes(&Q_, (uint8_t *)Q);
    peace_sm2_jacobian_point_add(&P_, &P_, &Q_);
    peace_sm2_jacobian_point_to_bytes(&P_, (uint8_t *)R);

    return 1;
}

/**
* @brief Subtracts one point from another.
* @param[out] R Result of the subtraction.
* @param[in] P Point to subtract from.
* @param[in] Q Point to subtract.
* @pre P and Q must be valid points.
* @post R will contain the result of the subtraction.
* @return 1 if successful, E_ERROR otherwise.
*/
int peace_sm2_point_sub(peace_sm2_POINT *R, const peace_sm2_POINT *P, const peace_sm2_POINT *Q)
{
    peace_sm2_JACOBIAN_POINT P_;
    peace_sm2_JACOBIAN_POINT Q_;

    peace_sm2_jacobian_point_from_bytes(&P_, (uint8_t *)P);
    peace_sm2_jacobian_point_from_bytes(&Q_, (uint8_t *)Q);
    peace_sm2_jacobian_point_sub(&P_, &P_, &Q_);
    peace_sm2_jacobian_point_to_bytes(&P_, (uint8_t *)R);

    return 1;
}

/**
* @brief Negates a point.
* @param[out] R Result of the negation.
* @param[in] P Point to negate.
* @pre P must be a valid point.
* @post R will contain the negated point.
* @return 1 if successful, E_ERROR otherwise.
*/
int peace_sm2_point_neg(peace_sm2_POINT *R, const peace_sm2_POINT *P)
{
    peace_sm2_JACOBIAN_POINT P_;

    peace_sm2_jacobian_point_from_bytes(&P_, (uint8_t *)P);
    peace_sm2_jacobian_point_neg(&P_, &P_);
    peace_sm2_jacobian_point_to_bytes(&P_, (uint8_t *)R);

    return 1;
}

/**
* @brief Doubles a point.
* @param[out] R Result of the doubling.
* @param[in] P Point to double.
* @pre P must be a valid point.
* @post R will contain the doubled point.
* @return 1 if successful, E_ERROR otherwise.
*/
int peace_sm2_point_dbl(peace_sm2_POINT *R, const peace_sm2_POINT *P)
{
    peace_sm2_JACOBIAN_POINT P_;

    peace_sm2_jacobian_point_from_bytes(&P_, (uint8_t *)P);
    peace_sm2_jacobian_point_dbl(&P_, &P_);
    peace_sm2_jacobian_point_to_bytes(&P_, (uint8_t *)R);

    return 1;
}

/**
 * @brief Multiplies a point on the elliptic curve by a scalar.
 *
 * This function multiplies the point P on the elliptic curve by the scalar k
 * and stores the result in R.
 *
 * @param[out] R Pointer to the resulting point after multiplication.
 * @param[in] k Pointer to the 32-byte scalar.
 * @param[in] P Pointer to the point to be multiplied.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_point_mul(peace_sm2_POINT *R, const uint8_t k[32], const peace_sm2_POINT *P)
{
	peace_sm2_BN _k;
	peace_sm2_JACOBIAN_POINT _P;

	// Convert the scalar k from bytes to a BN structure
	peace_sm2_bn_from_bytes(_k, k);
	// Convert the point P from bytes to a Jacobian point structure
	peace_sm2_jacobian_point_from_bytes(&_P, (uint8_t *)P);
	// Multiply the Jacobian point by the scalar k
	peace_sm2_jacobian_point_mul(&_P, _k, &_P);
	// Convert the resulting Jacobian point back to bytes
	peace_sm2_jacobian_point_to_bytes(&_P, (uint8_t *)R);

	// Clear the scalar k from memory
	memset(_k, 0, peace_sm2_sizeof_bn);
	return 1;
}

/**
 * @brief Multiplies the generator point on the elliptic curve by a scalar.
 *
 * This function multiplies the generator point G on the elliptic curve by the scalar k
 * and stores the result in R.
 *
 * @param[out] R Pointer to the resulting point after multiplication.
 * @param[in] k Pointer to the 32-byte scalar.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_point_mul_generator(peace_sm2_POINT *R, const uint8_t k[32])
{
	peace_sm2_BN _k;
	peace_sm2_JACOBIAN_POINT _R;

	// Convert the scalar k from bytes to a BN structure
	peace_sm2_bn_from_bytes(_k, k);
	// Multiply the generator point by the scalar k
	peace_sm2_jacobian_point_mul_generator(&_R, _k);
	// Convert the resulting Jacobian point back to bytes
	peace_sm2_jacobian_point_to_bytes(&_R, (uint8_t *)R);

	// Clear the scalar k from memory
	memset(_k, 0, peace_sm2_sizeof_bn);
	return 1;
}

/**
 * @brief Multiplies a point on the elliptic curve by a scalar and adds another point multiplied by another scalar.
 *
 * This function multiplies the point P on the elliptic curve by the scalar k,
 * multiplies the generator point G by the scalar s, and adds the two results,
 * storing the final result in R.
 *
 * @param[out] R Pointer to the resulting point after multiplication and addition.
 * @param[in] k Pointer to the 32-byte scalar for point P.
 * @param[in] P Pointer to the point to be multiplied by k.
 * @param[in] s Pointer to the 32-byte scalar for the generator point G.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_point_mul_sum(peace_sm2_POINT *R, const uint8_t k[32], const peace_sm2_POINT *P, const uint8_t s[32])
{
	peace_sm2_BN _k;
	peace_sm2_JACOBIAN_POINT _P;
	peace_sm2_BN _s;

	// Convert the scalar k from bytes to a BN structure
	peace_sm2_bn_from_bytes(_k, k);
	// Convert the point P from bytes to a Jacobian point structure
	peace_sm2_jacobian_point_from_bytes(&_P, (uint8_t *)P);
	// Convert the scalar s from bytes to a BN structure
	peace_sm2_bn_from_bytes(_s, s);
	// Multiply the Jacobian point P by k and add the generator point G multiplied by s
	peace_sm2_jacobian_point_mul_sum(&_P, _k, &_P, _s);
	// Convert the resulting Jacobian point back to bytes
	peace_sm2_jacobian_point_to_bytes(&_P, (uint8_t *)R);

	// Clear the scalars k and s from memory
	memset(_k, 0, peace_sm2_sizeof_bn);
	memset(_s, 0, peace_sm2_sizeof_bn);
	return 1;
}

/**
 * @brief Converts a point on the elliptic curve to compressed octets.
 *
 * This function converts the point P to a compressed octet representation
 * and stores the result in out.
 *
 * @param[in] P Pointer to the point to be converted.
 * @param[out] out Pointer to the buffer where the compressed octets will be stored.
 */
void peace_sm2_point_to_compressed_octets(const peace_sm2_POINT *P, uint8_t out[33])
{
	// Set the first byte to 0x02 or 0x03 based on the least significant bit of the y-coordinate
	*out++ = (P->y[31] & 0x01) ? 0x03 : 0x02;
	// Copy the x-coordinate to the output buffer
	memcpy(out, P->x, 32);
}

/**
 * @brief Converts a point on the elliptic curve to uncompressed octets.
 *
 * This function converts the point P to an uncompressed octet representation
 * and stores the result in out.
 *
 * @param[in] P Pointer to the point to be converted.
 * @param[out] out Pointer to the buffer where the uncompressed octets will be stored.
 */
void peace_sm2_point_to_uncompressed_octets(const peace_sm2_POINT *P, uint8_t out[65])
{
	// Set the first byte to 0x04 (uncompressed format)
	*out++ = 0x04;
	// Copy the point P to the output buffer
	memcpy(out, P, 64);
}

/**
 * @brief Converts octets to a point on the elliptic curve.
 *
 * This function converts the octets in the input buffer to a point on the elliptic curve
 * and stores the result in P.
 *
 * @param[out] P Pointer to the resulting point.
 * @param[in] in Pointer to the input buffer containing the octets.
 * @param[in] inlen Length of the input buffer.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_point_from_octets(peace_sm2_POINT *P, const uint8_t *in, size_t inlen)
{
	if ((*in == 0x02 || *in == 0x03) && inlen == 33)
	{
		// Convert from compressed octets
		if (peace_sm2_point_from_x(P, in + 1, *in) != 1)
		{
			return E_ERROR;
		}
	}
	else if (*in == 0x04 && inlen == 65)
	{
		// Convert from uncompressed octets
		if (peace_sm2_point_from_xy(P, in + 1, in + 33) != 1)
		{
			return E_ERROR;
		}
	}
	else
	{
		// Invalid input format
		return E_ERROR;
	}
	return 1;
}

/**
 * @brief Converts DER-encoded octets to a point on the elliptic curve.
 *
 * This function converts the DER-encoded octets in the input buffer to a point on the elliptic curve
 * and stores the result in P.
 *
 * @param[out] P Pointer to the resulting point.
 * @param[in,out] in Pointer to the input buffer containing the DER-encoded octets.
 * @param[in,out] inlen Length of the input buffer.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_point_from_der(peace_sm2_POINT *P, const uint8_t **in, size_t *inlen)
{
	int ret;
	const uint8_t *d;
	size_t dlen;

	// Extract the octet string from the DER-encoded input
	if ((ret = peace_asn1_octet_string_from_der(&d, &dlen, in, inlen)) != 1)
	{
		if (ret < 0)
			return ret;
	}
	// Check if the length of the octet string is 65 bytes
	if (dlen != 65)
	{
		return E_ERROR;
	}
	// Convert the octet string to a point on the elliptic curve
	if (peace_sm2_point_from_octets(P, d, dlen) != 1)
	{
		return E_ERROR;
	}
	return 1;
}

/**
 * @brief Converts a hash to a point on the elliptic curve.
 *
 * This function converts the hash of the input data to a point on the elliptic curve
 * and stores the result in R.
 *
 * @param[out] R Pointer to the resulting point.
 * @param[in] data Pointer to the input data to be hashed.
 * @param[in] datalen Length of the input data.
 * @return 1 on success, E_ERROR on failure.
 */
int peace_sm2_point_from_hash(peace_sm2_POINT *R, const uint8_t *data, size_t datalen)
{
	peace_sm2_BN u;
	peace_sm2_Fp x;
	peace_sm2_Fp y;
	peace_sm2_Fp s;
	peace_sm2_Fp s_;
	uint8_t dgst[32];

	// Calculate u = (P + 1) / 4
	peace_sm2_bn_add(u, peace_sm2_P, peace_sm2_ONE);
	peace_sm2_bn_rshift(u, u, 2);

	do
	{
		// Compute the hash of the input data
		peace_sm3_digest(data, datalen, dgst);

		// Convert the hash to a BN structure
		peace_sm2_bn_from_bytes(x, dgst);
		// Ensure x is less than P
		if (peace_sm2_bn_cmp(x, peace_sm2_P) >= 0)
		{
			peace_sm2_bn_sub(x, x, peace_sm2_P);
		}

		// Calculate s = (x^3 + ax + b) mod P
		peace_sm2_fp_sqr(s, x);
		peace_sm2_fp_sub(s, s, peace_sm2_THREE);
		peace_sm2_fp_mul(s, s, x);
		peace_sm2_fp_add(s, s, peace_sm2_B);

		// Calculate y = s^((P + 1) / 4) mod P
		peace_sm2_fp_exp(y, s, u);
		// Calculate s_ = y^2 mod P
		peace_sm2_fp_sqr(s_, y);

		// Update data and datalen for the next iteration
		data = dgst;
		datalen = sizeof(dgst);

	} while (peace_sm2_bn_cmp(s, s_) != 0);

	// Convert the x and y coordinates to bytes and store in R
	peace_sm2_bn_to_bytes(x, R->x);
	peace_sm2_bn_to_bytes(y, R->y);
	return 1;
}
